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Abstract 

Localization of plastic strain induced by softening can be objectively described 
by a regularized plasticity model that postulates a dependence of the current yield 
stress on a nonlocal softening variable defined by a differential (gradient) expression. 
This paper presents analytical solutions of the one-dimensional localization problem 
under certain special nonuniform stress distributions. The one-dimensional problem 
can be interpreted as describing either a tensile bar with variable cross section, or a 
beam subjected to a nonuniform bending moment. Explicit as well as implicit gradient 
formulations are considered. The evolution of the plastic strain profile and the shape 
of the load-displacement diagram are investigated. It is shown that even if the local 
constitutive law exhibits softening right from the onset of yielding, the global load- 
displacement diagram has a hardening part. The interplay between the internal length 
scales characterizing the material and the geometry is discussed. 

1 Introduction 

1.1 Gradient plasticity as a localization limiter 

Stress-strain diagrams of quasibrittle materials typically exhibit softening, which is caused by 
propagation and coalescence of defects in the microstructure, and at the macroscopic scale is 
manifested by decreasing stress at increasing strain. Softening can naturally be incorporated 
into damage models but, alternatively, can be modeled within the framework of plasticity 
with yield limit degradation. 

In the context of the standard continuum approach, softening may lead to localization 
of plastic strain into an arbitrarily small volume and, consequently, to the pathological 
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sensitivity of the numerical results obtained by the finite element method to the size of 
elements used in the analysis. Objective, mesh-insensitive description of the localization 
phenomenon requires an enhancement of the governing equations by non-standard terms 
that act as localization limiters and regularize the solutions. One popular class of localization 
limiters is based on the incorporation of second or higher gradients of the softening variable 
into the softening law [2]. Many gradient plasticity formulations have been inspired by the 
pioneering work of Aifantis and coworkers jTJ , others by the implicit gradient approach used 
in damage mechanics |10j . 

Localization properties of gradient plasticity models have often been studied by investi- 
gating the bifurcation from a uniform state in the idealized one-dimensional setting; see [5] 
for a summary. The aim of this paper is to extend the localization analysis to cases with a 
non-uniform stress field, which arises for instance in a bar under uniaxial tension due to a 
variation of the sectional area. To make the problem amenable to an analytical solution, we 
restrict our interest to static equilibrium under vanishing body forces, using the small-strain 
assumptions. 

1.2 One-dimensional softening plasticity model 

In the one-dimensional setting, the standard elastoplastic model is based on the additive 
split of the total strain e into the elastic part e e and plastic part e p , with stress a linked to 
the elastic strain by Hooke's law 

a = Ee e = E(e - e p ) (1) 

The yield function that identifies elastic, plastic and inadmissible states can be defined aa3 

f(a, k) = \a\ - oy(/e) (2) 

where ay is the current yield stress, evolving as a function of an internal variable n, to be 
specified later. The hardening-softening law will be considered in the simplest linear form, 

oy («0 = 0"o + Hk (3) 

where cq is the initial yield stress and H is the plastic modulus. Positive values of H 
correspond to hardening and negative values to softening. Our attention will be focused on 
the latter case, in which H is referred to as the softening modulus. 
The evolution of plastic strain is formally described by the flow rule 

e p = A sgn cr (4) 

lr The fact that we assume a symmetric response in tension and in compression, and that the hardening 
or softening is taken as isotropic, does not impose any restrictions on generality, since we will investigate 
localization under monotonic loading and the sign of stress or of the plastic strain increment will never 
change. 
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and the loading-unloading conditions 

A>0, /(<7,«)<0, A/(<7,k) = (5) 

where A is the plastic multiplier and the dot over a symbol denotes differentiation with 
respect to time. The internal variable k is usually taken as the cumulative plastic strain and 
defined by the rate equation 

k = \i p \ (6) 

If we consider only tensile loading (with possible elastic unloading, but never with a reversal 
of plastic flow), then the plastic strain e p , cumulative plastic strain k and plastic multiplier 
A are all equal. We will use k as the primary symbol for (cumulative) plastic strain and 
rewrite equations (pQ) and (JSJ) as 

a = E(e- k) (7) 

k>0, /(<7,k)<0, k/(<t,k) = (8) 

Formally the same framework can be used for the description of an elastoplastic moment- 
curvature relation that characterizes the inelastic flexural response of a beam. Stress and 
strain are then replaced by bending moment and curvature, Young's modulus E by the 
sectional bending stiffness EI (where I is the sectional moment of inertia), the initial yield 
stress o"o by the initial plastic moment Mo, and the softening modulus if by a constant C that 
represents the derivative of the plastic moment with respect to the plastic part of curvature 
(in the post-yield range). Description of the moment- curvature relation by a bilinear diagram 
is certainly a rough approximation, but it can reflect the main features of inelastic bending 
and serve as a prototype model, for which analytical solutions exist. 

1.3 One-dimensional localization problem 

It is well known that, in the one- dimensional setting, softening immediately leads to local- 
ization of inelastic strain. If we consider a straight bar with perfectly uniform properties 
subjected to uniaxial tension (induced by applied displacement at one bar end), the response 
remains uniform in the elastic range and also during plastic yielding with a positive plastic 
modulus. For a negative (or vanishing) plastic modulus, uniqueness of the solution is lost 
right at the onset of softening (or of perfectly plastic yielding). Stress distribution along 
the bar must still remain uniform due to the static equilibrium conditions (in the absence 
of body forces), but a given stress level can be attained by softening with increasing plastic 
strain, or by elastic unloading with no plastic strain evolution. Which cross sections un- 
load and which exhibit softening remains completely undetermined, and there is no lower 
bound on the total length of the softening region(s). Therefore, infinitely many solutions 
exist, including solutions with plastic strain evolution localized into extremely small regions. 
Since the dissipation per unit volume is fixed and the volume of the softening material is 
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arbitrarily small, failure of the bar can occur at arbitrarily small dissipation. This theoret- 
ical deficiency of the standard softening model is at the origin of numerical problems with 
pathological sensitivity of finite element solutions to the size of the elements. Even if the 
nonuniqueness of the solution is removed by a slight perturbation of the perfect uniformity 
of the bar, the problem with localization of softening into arbitrarily small regions (in fact 
into the weakest cross section) still persists. 

Remedy is sought either in adjustment of the softening modulus proportionally to the el- 
ement size [T2l [3] , or in advanced regularization techniques that serve as localization limiters, 
i.e., prevent localization of plastic strain into arbitrarily small regions and enforce a certain 
minimum size of the plastic zone. Of course this minimum size is not directly prescribed, it 
is rather the outcome of the solution of the governing equations, which are enhanced by non- 
standard terms that contain (sometimes in a hidden form) at least one new model parameter 
with the dimension of length. Such spatial scale information is not contained in standard 
constitutive equations based on traditional continuum mechanics. For instance, all the ma- 
terial parameters used by the simple elastoplastic model in subsection II. 2\ i.e. E, a and H, 
have the dimension of stress and cannot be combined into a parameter with the dimension 
of length. This is generally true of all the "local" and "standard-order" continuum theories^] 
Roughly speaking, by "local" we mean that the constitutive response at each material point 
depends on the state variables at that point only (and possibly on their previous evolution), 
and by "standard-order" we mean that (i) the state of deformation at each material point 
is fully described by the first gradient of the displacement field while the higher gradients 
have no influence, (ii) no additional kinematic variables are introduced, and (iii) gradients 
of internal variables are not incorporated in the constitutive equations. 

Enrichment of the constitutive model by nonstandard terms that contain a length pa- 
rameter can be achieved either by nonlocal approaches (in the narrow sense), which consider 
the influence of a finite neighborhood of each material point, or by higher-order approaches, 
which consider the influence of higher (than usual) gradients of the displacement or of internal 
variables, or introduce additional variables characterizing the deformation state (e.g. micro- 
rotations that are not equal to the macro- rotations derived from the displacement field), or 
deal with gradients of internal variables. Approaches based on higher-order gradients are 
sometimes considered as nonlocal in the broad sense, and certain specific formulations of this 
kind are even equivalent to nonlocal approaches in the narrow sense (which are also called 
integral-type nonlocal approaches because they use weighted spatial averaging described by 
integral operators). 

This study focuses on softening plasticity models enhanced by gradients of internal vari- 
ables, which means, for the simple one- dimensional model considered here, models with 
spatial derivatives of the cumulative plastic strain k. Localization properties of a wide range 
of such models have been analyzed and compared in [8], but only for the highly idealized 



2 A characteristic length is present in fracture mechanics, but that is not a "pure" continuum theory, since 
it admits discontinuities in the displacement field. 
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Figure 1: Tensile bars with (a) smooth distribution of sectional area, (b) continuous but 
non-smooth distribution of sectional 

case of a perfectly uniform stress field. It has been shown that some formulations suffer by 
serious deficiencies and thus do not need to be considered in more detailed studies. For this 
reason, we will restrict our attention to two widely popular families of models, referred to as 
explicit and implicit. Explicit gradient models have their origins in the pioneering work of 
Aifantis [lj , while implicit gradient models found inspiration in the implicit gradient damage 
model [10] and for plasticity were first proposed by Geers and coworkers [3, EJ [6] . In each 
family, we will investigate the basic version and one of its modifications: 

1. Explicit gradient model 

(a) with second gradient of cumulative plastic strain, 

(b) with fourth gradient of cumulative plastic strain. 

2. Implicit gradient model 

(a) with boundary conditions at the physical boundary, 

(b) with "boundary" conditions at the boundary of the plastic zone. 

In contrast to [SJ, we will consider nonuniform stress distribution, which can be caused 
by variations of the sectional area in the problem of bar under uniaxial tension, but at the 
same time can reflect the nonuniform distribution of bending moments in a beam. From this 
point of view, two basic cases can be distinguished: 

1. Smooth variation of stress, approximated in the vicinity of the global maximum by a 
concave quadratic function. 

2. Continuous variation of stress with discontinuous spatial derivative at the global max- 
imum, approximated by a concave piecewise linear function. 

The least regular case of a stress field with discontinuities (due to jumps in the sectional 
area) is more intricate and will be addressed in a separate publication. 

The first case corresponds to an axially stretched bar of a dog-bone shape; see Fig. [Ih, 
The stress distribution that allows for an analytical solution is given by the quadratic function 
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Figure 2: Simply supported beam with (a) quadratic distribution of bending moments, (b) 
linear distribution of bending moments 

where the origin of the spatial coordinate x is placed at the weakest section with maximum 
stress a c = cr(0), and l g = ^/— 2er(0)/er"(0) is a parameter that sets the length scale of stress 
variation. Expression (Q can be considered as the Taylor expansion truncated after the 
quadratic term. It would be exact for the special case of a bar with cross-sectional area 
varying according to 

«*) = 7^ = 7^ (10) 

where A c is the area of the weakest section and l g needs to be larger than the distance of 
that section from the bar end. However, even for bars with a more general but smooth vari- 
ation of cross-sectional area, the quadratic stress distribution (Q is a good approximation 
of the actual one in the vicinity of the weakest section, where the plastic zone is expected 
to develop. Therefore, analytical solutions of this special case can be considered as repre- 
sentative of other cases with more general but still smooth stress variations. At the same 
time, if the mathematical problem is interpreted as describing a bending beam instead of 
an axially loaded bar, the quadratic distribution of bending moments exactly corresponds 
to the fundamental case of a simply supported beam of length 2l g subjected to a uniform 
transversal load; see Fig. [2k 

The second continuous variation of stress with discontinuous spatial derivative, 

corresponds to an axially stretched bar with a wide V-shaped notch; see Fig. [lb. The stress 
distribution to be used in the analytical solution is given by 

*(x)=a c (l-^j (11) 

Again, for this relation to be exact, the notch would need to have a specific shape, with the 
sectional area varying according to 

M x ) = -4 = r% (12) 

1 _ M Ig-\X\ 
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From a more general point of view, the simple piecewise linear stress distribution serves as 
a prototype of all distributions with a kink at the weakest section. In terms of the beam 
bending problem, it exactly corresponds to another fundamental case of a simply supported 
beam of length 2l g subjected to a concentrated transversal force at midspan; see Fig. [2b. 

For each type of stress distribution and each gradient plasticity model, the post-yield re- 
sponse will be analyzed in terms of the distribution of plastic strain, evolution of plastic zone 
size, and load-displacement diagram. This will provide insight into the interplay between 
the material length scale introduced by the gradient enhancement of the elastoplastic model 
and the geometric length scale related to the variation of sectional area or, for the bending 
problem, to the span of the beam. 

2 Explicit second-order gradient plasticity model 

The explicit gradient plasticity model directly incorporates the second spatial gradient of 
cumulative plastic strain into the softening law. In the case of linear softening and in the 
one-dimensional setting, equation is replaced by 



where I is a new material parameter with the dimension of length, and primes denote deriva- 
tives with respect to the spatial coordinate x. In a general multidimensional setting, the 
second spatial derivative would be replaced by the Laplace operator. 

In the elastic range, the plastic strain identically vanishes and the yield stress cry is at 
its initial level, &q. As long as the stress cr(x) is everywhere below a , the response must 
remain elastic. The onset of plastic yielding occurs when the stress in the weakest section, 
a c , reaches the initial yield stress, oo. For a standard (not enriched) softening model, plastic 
yielding would localize into this single section. The gradient enhancement is supposed to 
act as a localization limiter, which makes the plastic zone grow to a finite size L p . This 
process will now be studied analytically for the two basic types of stress distributions, (j5j) 
and (fTTj) . It can be expected that a contiguous plastic zone forms around the weakest section 
and the surrounding parts of the bar remain elastic. Due to symmetry, the plastic zone is 
assumed to be an interval I p = (— L p /2, L p /2), symmetric with respect to the origin. Of 
course, the boundaries between the elastic and plastic zones can move, but as long as the 
plastic zone does not shrink (i.e., L p grows monotonically) , the elastic zone is characterized 
by zero plastic strain. In the plastic zone, the yield function must vanish, which provides 
an equation for the determination of plastic strain. For the gradient-enhanced model, this 
equation has a differential character. 




(13) 
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2.1 Quadratic stress distribution 

At all points of the plastic zone, the current yield stress ay given by (|T3|) must be equal to 
the stress a. For the quadratic stress distribution this condition leads to the equation 

k{x) + 1 2 k"{x) = - 1 (14) 

H 

This is a nonhomogeneous second-order linear differential equation with constant coefficients, 
and its general solution can be constructed as the sum of a particular solution of the nonho- 
mogeneous equation and the general solution of the corresponding homogeneous equation. 
For the quadratic function on the right-hand side of ([141) . there exists a quadratic particular 
solution 

k(x) = Ai + A 2 x 2 (15) 

with constants Ai and A 2 easily identified by substituting ffTBI) into the left-hand side of (fT41) 
and comparing the constant terms and the quadratic terms on both sides. By adding the 
general solution of the homogeneous equation, which is a linear combination of cos(x/l) and 
sin(x//), the general solution of (fT4"I) is obtained in the form 



a c (l 2 q + 2l 2 -x 2 )-a lt 



x „ . x 



= --vg ; -__ > -»g + Ci CQS - + C2 sin - (16) 

Integration constants C% and C 2 need to be determined from conditions k = and k,' = 
imposed at the boundary of the plastic zone. These are sometimes considered as boundary 
conditions, but they are better justified by regularity requirements. The plastic strain iden- 
tically vanishes outside the plastic zone, and if it did not remain continuously different iable 
across the elasto-plastic boundary, its second derivative would have a singular character and 
the yield condition would be "strongly violated". Precise mathematical justification of these 
statements, based on a variational formulation of the problem, will be presented in a sepa- 
rate paper. For the present purpose it is sufficient to admit that the plastic strain and its 
spatial derivative must vanish at the boundary of the interval that corresponds to the plastic 
zone. The end points of this interval are not known in advance, so in total we would have 
four unknowns (two integration constants and two coordinates of the end points) and four 
conditions. 

However, due to symmetry of the problem, the plastic zone I p = (—L p /2, L p /2) is centered 
at the origin, and the plastic strain distribution is described by an even function. This latter 
condition implies that integration constant C 2 must vanish. So it remains to determine C\ 
and L p from conditions 

k{L p /2) = 0, k'{L p /2) = (17) 
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Substituting the general solution ([IS]) with Ci = into (IT7|) leads to two equations 

a c (/ g 2 + 2Z 2 - \L 2 p ) - a l 2 g + C 1 Hl 2 g cos ^ = (18) 

or c L p Z + d#/ 2 sin ^ = (19) 

that are linear in terms of C\ and nonlinear in terms of L p . Unknown C\ is easily eliminated 
and the resulting equation for unknown L p reads 

« = a c (P 9 + 21°- \LI) - a P g (20) 

Recall that parameter o"o is a fixed material property while parameter a c represents the 
current stress in the weakest section and is directly related to the axial force transmitted by 
the bar, F = A c a c . Plastic yielding starts when <j c = ctq, i.e., when F = F where F = A c ao 
is the limit elastic force. For a given axial force F, the corresponding size of the plastic zone 
L p could be obtained by solving nonlinear equation (1201) numerically, with a c set to F/A c . 
However, for some load levels there could be multiple solutions or none at all. It is much 
more convenient to revert the procedure and express the axial force in terms of the plastic 
zone size, because this can be done analytically: 



A c a l 2 

P g + 2/ 2 - \Ll - L p l cotan(L p /2/) 



= n , o;9 lr9 C ^—-77 ( 21 ) 



A better understanding of the role of individual parameters can be gained if the problem is 
described in terms of dimensionless quantities. The force or stress level is reflected by the 
load parameter 

*=£ = - (22) 

which is equal to one at the onset of yielding. The length variables l g and L p are taken 
relative to the material length parameter I. Therefore we introduce dimensionless variables 

K = yv Xg = l f (23) 

Variable X p is the ratio between one half of the plastic zone size and the characteristic 
material length. Factor 2 in the denominator simplifies the results, because the boundary of 
the plastic zone is at x = l\ p . Parameter X g is the ratio between the geometric and material 
characteristic length scales. In terms of the dimensionless variables, equation (1211) reads 
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dimensionless plastic zone size, X. 



Figure 3: Explicit second-order gradient plasticity model: Dependence between auxiliary 
parameters 5 2 and Si, respectively given by and (l4"2j) . and the dimensionless plastic zone 
size X p = L p /2l 




-3-2-10123 



dimensionless coordinate, x/l 



Figure 4: Explicit second-order gradient plasticity model, quadratic stress distribution: Evo- 
lution of plastic strain profile 
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Figure 5: Explicit second-order gradient plasticity model, quadratic stress distribution: Plas- 
tic part of load-displacement diagram — (a) close-up of the initial part, (b) complete diagram 



where 

S 2 = \ 2 p + 2A P cotanAp - 2 (25) 

is an auxiliary variable introduced for convenience — it permits isolating the effect of X g from 
the dependence between X p and 0. The relation between X p and S 2 is graphically represented 
by the solid curve in Fig. [31 Parameter 62 continuously grows from its initial value at 
X p = to its maximum value 52,max = 7r 2 /4 — 2 attained at X p = ir/2, and afterwards 
decreases and tends to minus infinity as X p approaches 7r from the left. The corresponding 
load parameter <fi (axial force transmitted by the bar normalized by its elastic limit value) 
continuously grows from <p — 1 at the elastic limit to its maximum value 

F A 2 
^ max F A 2 + 2-7tV4 1 ' 

and afterwards decreases to zero as A p approaches ir from the left. An important point is that, 
during this process, A p grows monotonically, and so the assumption of plastic zone expansion 
is verified and the solution is admissible. The loading process can be parameterized by A p 
ranging from to n. 

For each A p G [0, it) and the corresponding a c = <j)<Jo determined from (1241) . integration 
constant C\ can be computed from ( Tl9l) and the resulting expression 

d = ~^fh- (27) 
Hr sin A p 
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can be substituted into the general solution ([TBI) . This gives the function 



k(x) 



,(£ + 2/ 2 



2aJ 2 X i 
sin X p 



\p x 

cos — 



(28) 



describing the plastic strain distribution along the process zone. For the representation of 
the plastic strain profile, it is again convenient to use normalized variables. The plastic strain 
is already dimensionless, but it is useful to normalize it by the reference value Kf = —ao/H, 
which corresponds to the plastic strain at full softening to zero stress according to the 
local linear softening law. The spatial coordinate is naturally normalized by the material 
characteristic length I, which leads to the dimensionless coordinate £ = x/l. The distribution 
of normalized plastic strain, n n = n/nf, is then described by the function 



«n(0 



1 



/ 2\ 



\2 



sin A, 



cos£ + £ 2 -2 



(29) 



where <p depends on A p and A g according to 
of dimensionless plastic zone sizes X p , equation 
profiles. An example with A„ = 5 and A„ = 1, 2. 



For a fixed ratio X g = l g /l and a sequence 



(1291) provides a sequence of plastic strain 
2.5, 2.8, 3, 3.1 and vr is plotted in Fig. H 
It is clear that at each point the plastic strain evolves monotonically, which verifies the 
admissibility of our solution. 

Finally, integrating the plastic strain along the entire process zone, we obtain 



L p /2 



U, 



k(x) dx 



H 



-L p /2 



mi 



(30) 



Variable u p can be interpreted as the plastic part of the bar elongation. The total bar 
elongation is the sum of u p and the elastic elongation u e 



C e F where 



da; 



EA(x) 



(31) 



is the elastic compliance of the bar. In (1311) . £ denotes the one-dimensional domain repre- 
senting the entire bar (interval of length L). The elastic compliance depends on the total bar 
length L, and thus the overall load-displacement diagram in terms of F against u e + u p would 
also be affected by the bar length. Since the plastic part of elongation, u p , is independent of 
L (as long as the size of the plastic zone does not exceed the bar length, which could happen 
only for extremely short bars), the plastic part of the load-displacement diagram, expressed 
in terms of F against u p , depends only on parameters /, l g , H, ctq and A c . If we convert 
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it to the dimensionless representation, with u p divided by the positive constant Ikj and F 
divided by the elastic limit force Fq = <JoA c , we obtain 



For a fixed value of the dimensionless parameter X g = l g /l, equations (j2~4l and (1521 provide a 
parametric description of the plastic part of the load-displacement diagram, parameterized 
by X p . The geometric length l g is typically much larger than the material length I. In fact, 
for the above solution to be valid in the entire softening range, l g should be larger than one 
half of the maximum plastic zone size, irl, otherwise the boundary of the plastic zone would 
reach the end of the bar and the regularity conditions at x — ±L p /2 should be replaced by 
boundary conditions at x — ±L/2. Therefore, the diagrams in Fig. Oare plotted for X g = 20, 
10, 5 and 3.2. It is found that, in accordance with what was already seen in Fig. [5J, in the 
first stage of the yielding process the force transmitted by the bar continuously increases 
from the elastic limit value Fq to the maximum value F max = max -Fo and only later starts 
decreasing and eventually vanishes. Note that the limit X g — > oo corresponds to a bar with 
a uniform section, for which the plastic zone size jumps immediately to 2ttI and the global 
response does not exhibit hardening at all. 

Summarizing the results obtained in this subsection for the smooth nonuniform distribu- 
tion of stress and second-order explicit gradient plasticity model with linear softening, it can 
be concluded that plastic yielding initiates in the weakest cross section (with the maximum 
stress) and then progressively expands, first at increasing axial force. Even though the 
local constitutive law assumes softening right from the onset of plastic yielding, the contri- 
bution of the gradient enhancement makes the yield stress in the weakest section and its 
neighborhood grow, until a certain size of the process zone is reached. The maximum force 
is therefore higher than the limit elastic force, and the global load-displacement diagram 
exhibits hardening. From (12T)|) it follows that 



which means that the difference between the reciprocal value of the limit elastic load (i.e. 
the ultimate load for the standard model) and the reciprocal value of the actual ultimate 
load for the gradient model is proportional to the square of the ratio between the material 
length parameter I and the geometric length parameter l g . This observation provides some 
insight into the interplay between the length scales of the problem. For I approaching 
zero, the hardening effects fade away because the gradient regularization is suppressed and 
the model approaches a standard one. For l g approaching infinity, the hardening effect 
fades away because the bar approaches a prismatic one, with uniform stress distribution, for 
which the global response exhibits linear softening without any previous hardening and the 




(32) 




2 



(33) 
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uniform response bifurcates into solutions with plastic strain localized into an interval of size 
2tcI, the position of which is completely arbitrary. For the nonuniform bar with quadratic 
stress distribution, 2irl is the limit that the process zone size approaches as the axial force 
transmitted by the bar decreases to zero during the second stage of failure, characterized by 
global softening. 

The same discussion applies to problems of beam bending with quadratic distribution of 
bending moment, e.g., to the simply supported beam loaded by uniform transversal load; 
see Fig. (2h. The middle section transmitting the largest moment starts failing first, and for 
a standard model with softening incorporated into the moment-curvature relation, inelastic 
processes would localize into this single section and the beam would fail at zero dissipation. 
For the second-order explicit gradient model, with the current yield moment given by 



plastic yielding spreads to a finite segment of the beam, corresponding to an inelastic hinge. 
Note that the plastic part of curvature, k p , is now marked with a subscript p, because k 
is used for the total curvature, consisting of the elastic and plastic parts. The integral 
of k p along the plastic zone corresponds to the rotation of an equivalent inelastic hinge 
concentrated into a single cross section, same as u p in equation (13"0|) can be interpreted as 
the opening of an equivalent cohesive crack. The inelastic hinge is therefore first hardening 
and only later softening. The peak bending moment depends on the geometric length scale, 
in the case of beam bending on the beam span, which introduces a size effect: shorter beams 
appear to be stronger (in terms of the maximum bending moment that they can transmit). 

Solution of one half of the simply supported beam also applies to a cantilever, provided 
that a homogeneous Neumann boundary condition for the plastic curvature is enforced at 
the clamped support. The cantilever problem was analyzed in [I], with the conclusion that 
the explicit gradient model is inappropriate and cannot describe the evolution of the plastic 
hinge. However, this conclusion was based on the tacit assumption that the global response 
must be softening right from the onset of plastic yielding. The present analysis shows that 
a reasonable solution can be found, with a continuous evolution of all quantities and a very 
reasonable physical interpretation, once it is accepted that, under nonuniform stress, local 
softening does not necessarily imply global softening. It is quite natural that, for a gradient 
model that exhibits (at least in a weak sense) a nonlocal behavior, the strength of the struc- 
ture does not depend exclusively on the weakest section but is influenced by the distribution 
of local strength in a neighborhood of that section, the size of which is controlled by the 
material characteristic length. Here we mean strength in the broad sense, in the present 
example referring to the axial force or bending moment that can be transmitted by a cross 
section. We have considered a variation of this sectional strength due to a nonuniform sec- 
tional area, but qualitatively the same results would be obtained in the case of a nonuniform 
local material strength (initial yield stress a ). 




(34) 
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2.2 Piecewise linear stress distribution 

Now we proceed to the piecewise linear stress distribution from Fig. [lb, with a discontinuous 
spatial derivative at the weakest section. The analysis proceeds along similar lines as for 
the quadratic stress distribution in the previous subsection, and so we omit most of the 
detailed explanatory comments. Whenever possible, the solution will be presented in terms 
of dimensionless quantities. 

In the plastic zone, the stress given by ffTTj) must be equal to the yield stress given by 
( 1T3|) . which leads to the differential equation 

k(x) + / V(x) = - ^ (35) 



or, in dimensionless form, 



«;„(£) + <(£) = W + <M; X |£| (36) 



where n n = n/ftf, Kf = —a /H, = a c /cr , \ g = l g /l and, for simplicity, the derivative with 
respect to the dimensionless coordinate f = x/l is denoted by a prime, even though before 
the prime was used for the derivative with respect to x. The general solution of fl36l) is 

f 1 - (f> - (pX- 1 ^ + d cosf + C 2 sinf for - A p < f < 
I l-(f) + <j)\g~ £ + C 3 cos^ + C 4 sin^ for < f < A p 

Continuity at £ = implies that C1 — C3, and continuous differentiability at £ = combined 
with symmetry implies that 

C 2 = -C A = <p\- 1 (38) 



So the general solution fl37l) can be rewritten as 

K n {0 = l-4> + 0A; 1 (If I - sin If I) + d cosf (39) 

Integration constant Ci and the size of the plastic zone L p are determined from the conditions 
ftn(Ap) = and K' n {\) = 0. Elimination of the integration constant leads to a nonlinear 
equation 

6 (1 — cos A D ) , 
(0- l)A 9 + 0(sm Ap-Ap) 
that links the plastic zone size A p to the load parameter <f), with an influence of parameter 

The dependence on X g can be treated separately if we define an auxiliary parameter 

5 1 =(l- I) X 9 (41) 



A, 
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This parameter can be expressed from (|40l) as a function of the dimensionless plastic zone 
size: 

• 1 — cos Xp Xp , , 

01 = Xp — sin Xp H = Xp — tan — (42) 

tanAp 2 

The corresponding graph is represented by the dashed curve in Fig. [3j As X p grows, parameter 
Si first increases from its initial value 5i = at X p = to its maximum 5i, max = it/ 2 — 1 at 
X p = 7r/2, and then decreases and tends to minus infinity as X p approaches tt from the left. 
This means that, similar to the previous case of quadratic stress distribution, the axial force 
transmitted by the bar, 

I? — T? A. — F ° _ F ° A 9 (43) 



1 - tfiA" 1 X g - Xp + tan(A p /2) 



is first increasing, even after the onset of plastic yielding, and the global response is hardening 
up to the peak force 

Fr. A - F„ \ 

(44) 



FoXg F Q Xg 



A 9 — <5i jmax X g + 1 — 7r/2 

After that the force decreases to zero and the plastic zone keeps expanding up to its maximum 
size 2irl. 

Integration constant C\ can be expressed as 

CWA^tan^ (45) 
and substituted into which leads to the particular solution 



MO = 1 - + ( \i\ - sin \i\ + tan ^ cos£ ) (46) 



The normalized plastic displacement is then 



lK f J T 2A 



4>X P 



«n(0de = 2Ap 1-0 + ^ (47) 



V'/ 



The evolution of the plastic strain profile is shown in Fig. [6^ for A 9 = 5. Same as for 
the quadratic stress distribution, the plastic strain grows monotonically at each fixed spa- 
tial point, which verifies the admissibility of the solution. The plastic part of the load- 
displacement diagram is plotted in Fig. [7|for several values of X g . 

It is also instructive to look at the distribution of the plastic strain rate. In equation 
15]) . A j is a fixed parameter and £ is the spatial coordinate. Only the plastic zone size X p 
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Figure 7: Explicit second-order gradient plasticity model, piecewise linear stress distribution: 
Plastic part of load-displacement diagram — (a) close-up of the initial part, (b) complete 
diagram 
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and the load parameter <p evolve in time. Their rates are linked by the rate form of equations 
(|3I]) -(H2j). from which 

02c ° sAp A n (48) 



2X g cos 2 4f 



Differentiating ( 1461) with respect to time and using ( )48l) . we obtain 



«n(f) 



If I -shi |f | +tan^cosf- \ g 



> cos X p + X g cos f 



0A, 



2A 2 cos 2 f 



(49) 



It is possible to show that k n (X p ) = but k r n (X p ) ^ 0. This means that the plastic strain rate 
is continuous but not continuously differentiable (in space) at the boundary of the plastic 
zone; see Fig. [6b. Therefore, the conditions of vanishing spatial derivative of plastic strain 
and of vanishing spatial derivative of plastic strain rate are in general not equivalent. Only 
in the special case of a plastic zone of constant width, which occurs in a bar with perfectly 
uniform properties, one can freely choose between conditions «4(A P ) = and k' n {X p ) = 0. 
The role of moving boundaries of the plastic zone has been examined for instance in |llj . 



3 Explicit fourth-order gradient plasticity model 

Modified versions of the second-order model from the previous section can incorporate higher- 
order gradients [13 E]. A prototype model of this kind, using a fourth-order enrichment, 
postulates the softening law in the form 

a Y = a + H (« - 1 4 k iv ) (50) 

where superscript "/V" denotes the fourth spatial derivative. In a general multiaxial setting 
this would be the "Laplacean of the Laplacean" . 



3.1 Quadratic stress distribution 

For the quadratic stress distribution (jUJ), the distribution of plastic strain inside the plastic 
zone I p = (—Lp/2, L p /2) is governed by the fourth-order differential equation 

. ri/ cr c — crn — cr c x 2 // 2 

k(x) - ZV y (x) = — (51) 

H 

which can be obtained by combining the yield condition o = oy with the softening law ( 1501) . 
In the dimensionless format, we rewrite it as 



MO->4 y (O = i-0 + ^r (52) 
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and construct the general solution 

<hf 2 

«n(0 = 1 - <\> + -^r + Ci cos£ + C 2 sin£ + C 3 cosh£ + C 4 sinh£ (53) 
\ 

where, as usual, K n — k/k/ — —kH/(To, £ = x/l, = <j c /ijq and X g = l g /l. 

Integration constants C 2 and C4 must vanish because of symmetry. The remaining un- 
knowns are integration constants C\ and C3 and the size of the plastic zone L p = 21 X p , and 
they need to be determined from regularity conditions at the boundary of the plastic zone, 
i.e., at point £ = X p . Since we have only three unknowns, we cannot impose continuous 
differentiability up to the third order, because this would represent four independent condi- 
tions. One may think that the discrepancy is caused by the assumption of symmetry. This 
is not the case — if symmetry is not imposed, we have to determine four integration constants 
and two coordinates of the boundary of the plastic zone (left and right boundary), which 
makes a total of six unknowns, and again only three independent conditions can be satisfied 
at each boundary point. Consequently, continuity at the boundary of the plastic zone can 
be enforced only for n n itself and its first and second derivative. As long as the plastic zone 
expands monotonically, the plastic strain is identically equal to zero outside the plastic zone, 
and the continuity conditions read 

K n (\ p ) = 0, <(A P ) = 0, <(Ap) = (54) 

The third derivative can be discontinuous, which means that the third derivative from the 
right (equal to zero) and the third derivative from the left, denoted as K r ^(X~), can be 
different. However, the jump in the third derivative, given by = — n r ^(X p ) = 

—k'^(X~), must be nonnegativejf] and so k'"(X~) must be nonpositive. This condition is to 
be verified once the solution is known. 

After substitution of the general solution (|53|) . conditions ( 1541) lead to the set of equations 

C\ cos X p + C 3 cosh X p = —l + (f> — (j)XpXg 2 (55) 
—G\ sin X p + C 3 sinh X p = — 24>X p X~ g 2 (56) 
—C\ cos Xp + C3 cosh X p = —2(pX g 2 (57) 

Elimination of integration constants C\ and C3 reduces the problem to a single nonlinear 
equation for X p . Using an auxiliary parameter 

8 2 = (l - I) X] (58) 

3 If the jump of the third derivative of plastic strain is positive, the fourth derivative contains a singular 
component that has the character of a positive multiple of the Dirac distribution. Since the plastic modulus 
H is negative, the yield stress at the boundary of the plastic zone computed from (|50|) is "positively infinite" 
and the solution remains plastically admissible. A precise justification of these rather intuitive arguments 
can be based on the reformulation of the problem as a variational inequality. 
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Figure 8: Explicit fourth-order gradient model: Dependence between dimensionless load 
parameters 5i and 5 2 and dimensionless plastic zone size X p 



the resulting equation can be presented in the dimensionless form 

4A P - (2 + A 2 -£ 2 )tanhA p 
tanA p = ^ (59) 



In a similar spirit as in Section [27TI instead of solving for A p , we can invert the problem and 

express 

4A P + (A 2 -2)tanA p - (A 2 + 2)tanhA p 



tan Ap — tanh A p 



(60) 



from which the load parameter is easily calculated as 



A 2 

(61) 



A 2 - 5 2 



9 



The advantage is that, using this representation, the influence of parameters A p and X g is 
treated separately. 

According to f|60l) . variable <5 2 continuously grows from its initial value at A p = to its 
maximum value <52,max ~ 0.779 at A p ~ 2.365, and then decreases and tends to minus infinity 
as Ap approaches 3.9266; see Fig. [SJ Expressing now the integration constants C\ and C 3 
from fl56l) - fl57j) and substituting into the general solution fl53l . we construct the particular 
solution 



cos^ cosh^ 

MO = i— 7s — x — o — r~r+ 

2 cos Ap 2 cosh A p 



e _ , ( x , cose L _ 2 + coshe ' 

X 2 g V X 2 g J 2cosAp ^ V A 2 / 2coshA p 

(62) 
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Figure 9: Explicit fourth-order gradient model, quadratic stress distribution: Evolution of 
(a) plastic strain, (b) third derivative of plastic strain 



The evolution of the plastic strain profile, graphically presented in Fig. [9^ for X g 
again monotonic. The third derivative of plastic strain, 



4, is 



'(0 



sinh£ 



2 cos A„ 2 cosh A r 



sin£ 



+ • 



\ 2 



_sinj_ 

2 cos A, 



+ 1 



2 + A; 



sinhtf 
2 cosh A r 



(63) 



is discontinuous at the boundary of plastic zone, i.e., at points £ = ±A P . This is documented 
by the graphs in Fig. [9b, plotted for X g = 4. The jumps in the third derivative are always 
positive, because A p ) > and k'^(X p ) < for each admissible plastic zone size X p , and 
so the solution is plastically admissible. 

Integrating the normalized plastic strain along the plastic zone, we obtain the normalized 
plastic elongation 



Up 



*2 



8(tanA p — A p )(tanhA p — A p 
tanh A p — tan A p 



4A 3 



(64) 



The resulting plastic part of the dimensionless load-displacement diagram is plotted in Fig. [10] 
for several values of parameter X g . 
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Figure 10: Explicit fourth-order gradient model, quadratic stress distribution: Plastic part 
of load-displacement diagram for different values of X g - - (a) early stages, (b) complete 
diagram 

3.2 Piecewise linear stress distribution 

For the piecewise linear stress distribution ffTTj) and the explicit fourth-order model, the yield 
condition leads to a differential equation that can be written in the dimensionless form as 



aO-C(£) 



+ 



A, 



(65) 



and has the general solution 

1 - 4> - <M~ X £ + Ci cos£ + C 2 sin£ + C 3 cosh£ + C A sinh£ for - X p < £ < 
1 - 4> + (f>Xj^ + C 5 cos£ + C 6 sin£ + C 7 cosh£ + C 8 sinh£ for < f < A p 



«n(f) 



(66) 

If symmetry is ignored, eight integration constants plus two coordinates of the boundary 
of the plastic zone can be determined from four continuity conditions at £ = and three 
continuity conditions at each of the boundary points. By exploiting symmetry, we can reduce 
the number of unknowns without affecting the solution. Symmetry conditions imply that 



W — ^5) ^2 — ^3 — ( - y 4 



-a 



(67) 



and they lead to an automatic satisfaction of continuity of K n and k'^ at £ = 0. We still need 
to enforce continuity of n' n and k'£ at £ = 0. From these conditions combined with symmetry 
we obtain 



Cfi — Cs 



A 

2A, 



(68) 
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The remaining three unknowns, C5, C7 and A p , need to be determined from the conditions of 
twofold continuous differentiability at £ = X p . Under the assumption of monotonic expansion 
of the plastic zone, these condition can be written in the same form ( 1541) as for the previous 
case with quadratic stress distribution. After substitution from ( 1661) and (IBS!) into ( 1541) . we 
obtain a set of three equations 

C 5 cos X p + Cj cosh A p = — 1 + -^-(sin A p + sinh A p — 2A P ) (69) 

2A 9 

—C5 sin Ap + C7 sinh Ap = ttt - (cos A p + cosh A p — 2) (70) 

2A g 

— C5 cos Ap + C7 cosh Ap = -^-(— sin Ap + sinh Ap) (71) 

2A 9 

and finally, after elimination of C5 and C7, we express the load parameter in the form 

A„ 



Xg~5l 

where the auxiliary parameter 



(72) 



5 1 = X + 2 ~ cos ^p-cosh 1 X P 
p tan A p — tanh A p 

depends only on parameter X p . The relation between £1 and X p is shown graphically by the 
dashed curve in Fig. El 

After evaluation of the integration constants and their substitution into (1661) . we obtain 
the particular solution 



cos £ cosli c" 

2|f | - 25i + (5 1 - Ap + sin Ap) + (5 1 -X p + sinh X p ) — — ^ - sin |£| - sinh |£| 

cos Ap cosh Ap 

(74) 

The plots in Fig. [TTh show a typical evolution of plastic strain, obtained for X g = 4. The plas- 
tic zone again starts from the weakest section and expands monotonically, first at increasing 
and later at decreasing axial force. The third derivative of plastic strain, 



sin sinli £ 

(5 1 - X p + sinA„) — + (5 1 — X p + sinhA„) — — — + sgn£(cos£ - cosh£) 

cosAp coshAp 



(75) 

is discontinuous at the boundary of plastic zone, i.e., at points £ = ±A P . This is documented 
by the graphs in Fig. [TTb . plotted for X g = 4. Same as in the case of quadratic stress dis- 
tribution, the jumps in the third derivative are always positive, and the plastic admissibility 
of the solution is verified. 
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Integrating the plastic strain (1741) . we obtain the plastic elongation, which can be pre- 
sented in the dimensionless form as 

— ^- = — [Ap — 2<5iA p + (Si — A p )(tan A p + tanh X p ) + cos -1 A p — cosh -1 A J (76) 

iKf Xg 

The plastic part of the load-displacement diagram is depicted in Fig. [12] for several values of 
parameter A g . 



4 Implicit gradient plasticity model 

The implicit gradient approach incorporates the influence of the nonlocal cumulative plastic 
strain R, defined as the solution of the differential equation 

R(x) - l 2 R"(x) = k{x) (77) 

usually with homogeneous Neumann boundary conditions 

R' = on dC (78) 

where £ denotes the interval representing the bar and dC is its boundary, consisting of two 
points. In multiple dimensions, the second spatial derivative is replaced by the Laplacean. 
The reason why R is called nonlocal is that its value at a given point x depends on the values 
of the "local" cumulative plastic strain k at all points of the body. It can even be shown 
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Figure 12: Explicit fourth-order gradient model, piecewise linear stress distribution: Plastic 
part of load-displacement diagram for different values of \ g — (a) early stages, (b) complete 
diagram 



that R corresponds to the weighted spatial average of n with a special choice of the weight 
function, set equal to the Green function of boundary value problem (1771) - (1781) . 

The simplest formulation of an implicit gradient plasticity model could be based on the 
replacement of the local cumulative plastic strain k in the softening law by its nonlocal 
counterpart R. This formulation enforces nonzero dissipation but does not prevent localiza- 
tion of plastic strain into a set of zero measure (into a single cross section). In analogy to 
integral-type nonlocal plasticity [HI [15] . a finite size of the process zone can be obtained 
with k replaced by an "overnonlocal" variable R, defined as the linear combination 

R = mR + (1 — vo)k (79) 

where parameter m is larger than 1. In such a case, the linear softening law is written as 

cr y = a + H[mR + (1 - m)/c] (80) 

In the plastic zone, the yield condition a = ay leads to the differential equation 

a(x) - a 



mK(x) + (1 — 7ti)k(x) 



H 



(81) 



with two unknown functions, R and k, which are linked by f[7T|) . Instead of solving a set of 
two differential equations, it is convenient to substitute (1771) into (18T1) and eliminate k. This 
leads to the equation 

(82) 



R(x) + (m - l)l 2 R"(x) 



a(x) 



H 
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with only one unknown function, R. 

In terms of the normalized functions K n (£) = k,(1^)/k,j and R n (Q — with 
Kf — —(To/H, equations (1771) and ( J82l) are rewritten as 

K n (x) - Kn( X ) = K n(x) (83) 

Ut) + V?<it) = (84) 

where the prime denotes differentiation with respect to the dimensionless coordinate £ = x/l, 
and /i = \/m — 1 is introduced for convenience. To get the simplest possible form of equa- 
tions, symbols m and /i will be used simultaneously, but they refer to a single independent 
material parameter, and m can always be replaced by 1 + /i 2 or \x by y/m — 1. 

Equation (1841) is valid only in the plastic zone, characterized by a positive value of the 
local variable K n . Outside the plastic zone, we have n n = and R n is governed by the 
homogeneous version of equation fl83|) . i.e., by 

M0-<(0 = o (85) 
4.1 Quadratic stress distribution 

For the quadratic stress distribution given by equation (1541) has the specific form 

MO+/i 2 <(£) = W+^V ( 86 ) 
where <fr — &c/ °"o is the load parameter. The general solution in terms of the nonlocal variable 

R n (0 = 1 - (f) ( 1 + 2/i ' 7 H + C lC os^ + C 2 sm- (87) 



substituted into (1831) gives the local variable 

«n(0 = MO -<(£)= 1 - (l + ^7^) + ^ (Ci cos | + C 2 sin 



By virtue of symmetry we get C 2 = 0. 

Equations fl86l) - fl88|) are valid only in the plastic zone I p characterized by n n > 0. Outside 
the plastic zone, the local variable k„ vanishes and the nonlocal variable R n is governed by 
the homogeneous differential equation ( |85l) with general solution 

Rn(0 = C 3 et + C 4 e-Z (89) 

At the physical boundary, the homogeneous Neumann boundary condition R' n = is usually 
imposed. For a finite bar, the solution then depends on the bar length as an additional 
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Figure 13: Implicit gradient model, quadratic stress distribution: Dependence between aux- 
iliary parameter <5 2 and dimensionless plastic zone size X p 

parameter. For simplicity, we assume that the bar is much longer than the plastic zone. The 
boundary condition is then imposed "at infinity", which means that it is replaced by the 
requirement that the solution must remain bounded. In that case, in the "right" part of the 
elastic zone integration constant C3 must vanish. Integration constants C\ and C4 and the 
dimensionless plastic zone size X p can be calculated from the condition n n (X p ) = and from 
the conditions of continuous differentiability of R n at £ = X p . 

After elimination of C\ and C4 we end up with a single equation that links the plastic 
zone size X p and the load parameter 0, of course with an influence of parameters X g and m. 
It turns out that, similar to (124]) and (I6T1) . the load parameter can be expressed in the form 



\ 2 



(90) 



but the auxiliary parameter 62 is instead of (1251) or fl60|) given by 



X v — utan — 

5 2 = X 2 p + 2m-^— ^ (91) 

1 + u tan — 

The dependence of 62 on X p is graphically illustrated in Fig. [13] for several values of parameter 
m. Once again, the advantage of using the auxiliary parameter is that parameter X g appears 
only in ( 190]) and parameters X p and m or /i only in ( 19T]) . The integration constant is then 
evaluated as 



20 (1 + X p )fi 



2 



C 1 = T2 A A ( 92 ) 

\ cos-f +//sm-f 



4 IMPLICIT GRADIENT PLASTICITY MODEL 



28 



(a) (b) 




-4 -2 2 4 -4 -2 2 4 



dimensionless coordinate, x/l dimensionless coordinate, x/l 



Figure 14: Implicit gradient model, quadratic stress distribution: Evolution of normalized 
plastic strain profiles for m = 4 and X g = 5 — (a) local plastic strain, (b) nonlocal plastic 
strain 



and the distribution of plastic strain in the process zone 



«n(0 



\2 



X p - 



2m(l + \ p ) 



1 + n tan 



1 - 




(93) 



is obtained by substituting f[9"2"j) and C 2 = into 

An example of the evolution of normalized local and nonlocal plastic strains for m = 4 
and X g = 5 is presented in Fig. [TH Note that the local plastic strain is continuous but not 
continuously differentiable at the boundary of the plastic zone. In this aspect, the implicit 
gradient model differs from the explicit one. Outside the plastic zone, the local plastic strain 
vanishes but the nonlocal plastic strain does not. The evolution of local plastic strain is 
monotonic and the plastic zone expands up to its maximum length 



2A 



p,max 



2/i 



7T 



arctan — 



(94) 



which depends on parameter \x and tends to zero as fi approaches zero, i.e., as m = 1 + fi 2 
approaches 1 from the right; see Fig. [15j This confirms that the simple model with m — 1, 
i.e., with softening driven by the nonlocal cumulative plastic strain R, does not act as a 
genuine localization limiter. 

Integration of the normalized plastic strain leads to the normalized plastic elongation 



Mr, 



iKf 



\2 
X 9 



2m(l + X p ) 
1 + ri tan — 



/itan — — A p I — 



2X l 



(95) 
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Figure 15: Implicit gradient model: Dependence of the maximum dimensionless half-size of 
plastic zone, A Pjmax , on model parameter m 



The plastic part of the load-displacement diagram is shown in Fig. [16] for different values of 
parameter X g . For smaller values of X g , the peak load is higher and the dissipated energy as 
well. For large values of A 9 , the diagram is close to a linear one. 



4.2 Piecewise linear stress distribution 

For piecewise linear stress distribution pip , equation (|84p has the specific form 



«n(0+A:(0 = W+^ ( 96 ) 

For £ e [0, X p ], the general solution reads 

k„(0 = l-0 + ^ + C lC os- + C 2 sin^ (97) 

Xg H [L 

and from symmetry and continuous differentiability we get the condition R' n (0) = 0, which 
implies that 

C 2 = -H* (98) 

At the boundary of the plastic zone, we impose the condition of vanishing local plastic strain, 
^n(Ap) = 0, which can be rewritten in terms of the nonlocal plastic strain as K n (X p ) = R'^(X P ). 
Continuous differentiability of the nonlocal plastic strain combined with boundedness in the 
semiinfinite elastic zone leads to the condition R n (X p ) = —R' n (X p ). In this way, two equations 
for two unknowns C\ and X p are constructed. Eliminating the integration constant C±, we 
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Figure 16: Implicit gradient model, quadratic stress distribution: Plastic part of load- 
displacement diagram for m = 4 and different values of \ g — (a) early stages, (b) complete 
diagram 



arrive at one single equation that links the load parameter 
After some algebra, the load parameter can be expressed as 

A„ 



and the plastic zone size X p . 
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where 



Si = A p 



2m sin 
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is an auxiliary parameter that depends only on A p and m but not on \ g . The relation 
between 8\ and X p described by f llOOp is graphically presented in Fig. [17] for several values of 
parameter m. The plastic zone expands to its maximum size, given again by formula (194")) . 
which was graphically presented in Fig. [T5l 

For each admissible size of plastic zone A p , integration constant C% is readily expressed, 
and its substitution into (!97j) provides the particular solution. To make the resulting ex- 
pression valid in the entire plastic zone, we replace £ by its absolute value, accounting for 
symmetry. The final formula for the normalized nonlocal plastic strain is then 
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Figure 17: Implicit gradient model, piecewise linear stress distribution: Dependence between 
auxiliary parameter 6i and dimensionless plastic zone size A p 



and the corresponding normalized local plastic strain is given by 

MO = *n(£) -<(£) = 
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An example of evolution of the plastic strain profile is shown in Fig. [T8j Note the discontinuity 
in the derivative of the local plastic strain at £ = 0, which is clearly manifested at early stages 
of evolution. 

Integrating f 1 1 2 j) over the plastic zone, we obtain the dimensionless plastic elongation 
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Ap(A p — 25i 
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The corresponding load-displacement diagram is plotted in Fig. [19l 



(103) 



5 Implicit gradient plasticity with modified boundary 
conditions 

The usual formulation of implicit gradient plasticity, developed by Geers and coworkers 
[3, E] and considered in the previous section, imposes the homogeneous Neumann boundary 
condition at the physical boundary of the body of interest. In a recent study focusing 
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Figure 18: Implicit gradient model, piecewise linear stress distribution: Evolution of nor- 
malized plastic strain profiles for m = 4 and \ g = 5 — (a) local plastic strain, (b) nonlocal 
plastic strain 
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Figure 19: Implicit gradient model, piecewise linear stress distribution:Plastic part of load- 
displacement diagram for m = 4 and different values of X g — (a) early stages, (b) complete 
diagram 
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on applications to beam bending, Challamel [I] proposed to impose that condition on the 
boundary of the plastic zone. For comparison, we will present the solution of the one- 
dimensional localization problem based on this modification. The only difference compared 
to the analysis performed in the previous section is that the requirement of boundedness 
of the nonlocal plastic strain is replaced by the conditions K' n (X p ) = and R' n (—X p ) = 0. 
Since the analysis proceeds along the same line as in the previous section, we omit all the 
intermediate steps and directly proceed to the results and their discussion. 

5.1 Quadratic stress distribution 

The auxiliary parameter 5 2 related to the load parameter <fi by (190]) is expressed as 



5 2 = X 2 v + 2m f^cotan^-1 
and the distribution of nonlocal and local plastic strain is given by 



(104) 
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The plastic zone expands up to its maximum size 
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and the normalized plastic elongation is 
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The dependence of the auxiliary parameter 5 2 o n the plastic zone size is graphically 
shown in Fig. [20] for several values of model parameter m. Interestingly, if m is in the range 



between 1 and 3, parameter 5 monotonically decreases from its initial value at A p = and 
never becomes positive. This means that the corresponding load parameter (p given by f]9"U]) 
monotonically decreases from its initial value 1, and the global response does not exhibit 
any hardening at all; see Fig. [21] In this aspect, the current model differs from all the other 
models covered by the present study. Global hardening arises only for m > 3. 
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Figure 20: Modified implicit gradient model, quadratic stress distribution: Dependence 
between auxiliary parameter 5 2 and dimensionless plastic zone size \ p — (a) early stages, 
(b) global picture 
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Figure 21: Modified implicit gradient model, quadratic stress distribution: Plastic part of 
load-displacement diagram for different values of parameters m and \ g — (a) early stages, 
(b) complete diagram 
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Figure 22: Modified implicit gradient model, quadratic stress distribution: Evolution of 
normalized plastic strain profile for m = 2 and A s = 5 - - (a) local plastic strain, (b) 
nonlocal plastic strain 

The evolution of the local and nonlocal plastic strain profiles is depicted in Fig. [22] for 
the specific parameter values m = 2 and \ g = 5. Similar to all previously discussed models, 
the plastic zone expands and the local plastic strain grows monotonically. However, a deeper 
examination of the solution reveals that the current model exhibits a pathological behavior, 
which becomes apparent if we plot the distribution of the actual stress and the current yield 
stress (Fig. [23]). In the plastic zone I p that extends from £ = — X p to £ = X p , the actual 
stress a is equal to the yield stress cy, which is quite natural because the yield condition 
a = ay is in fact the differential equation from which the solution has been calculated. 
In the elastic zone, the solution should comply with the condition of plastic admissibility, 
requiring that the actual stress must not exceed the yield stress. As seen in Fig. [231 the 
inequality a < ay is satisfied in the proximity of the plastic zone but violated farther away. 
The problem is caused by the modified "boundary" condition, which enforces a vanishing 
derivative of nonlocal plastic strain at the boundary of the plastic zone. The nonlocal plastic 
strain has then local minima at the points separating the plastic zone from the elastic ones 
(Fig. 122b). and its value in the elastic zones, calculated from ([85]) . increases with increasing 
distance from the plastic zone (note the hyperbolic cosine function in the second line of 
( 11051) ). The local plastic strain in the elastic zones vanishes and the nonlocal plastic strain 
becomes large, which results into a dramatic reduction of the yield stress. So even though 
the actual stress in the elastic zone is relatively low, the plastic admissibility condition is 
violated and additional plastic zones would be formed. This is of course a non-physical, 
unacceptable artefact, and thus the model with the modified boundary condition cannot be 
considered as a viable alternative to the standard implicit gradient formulation. 
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Figure 23: Modified implicit gradient model, quadratic stress distribution: Distribution of 
normalized actual stress cr/a and yield stress cr Y /(r at plastic zone size (a) X p — 1, (b) 



X p — 2 



5.2 Piecewise linear stress distribution 

For completeness, the results obtained with the modified implicit gradient model for a piece- 
wise linear distribution of stress are presented: 

<Ji = A p --tan^ (109) 
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Since the model exhibits a pathological behavior, similar to the preceding subsection, the 
graphical presentation and discussion of the results can be omitted. 
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6 Summary and conclusions 

In this paper, analytical solutions have been derived for localization of plastic strain described 
by gradient plasticity models in one spatial dimension. In previous studies, this problem was 
usually investigated for the simplest case of a bar with uniform properties, and localization 
was treated as a bifurcation from the uniform solution. Here we have considered a nonuniform 
distribution of stress along the bar due to a variable cross-sectional area. To keep the problem 
tractable, the stress distribution has been assumed to be quadratic or piecewise linear. The 
first case represents a smooth variation of the sectional area and the second case represents 
a non-smooth but still differentiable variation. Alternatively, the resulting mathematical 
problem can be interpreted as the description of a bending beam subjected to a uniform or 
concentrated lateral load. The function specifying the stress or bending moment distribution 
contains a parameter with the dimension of length, which reflects the length scale of the 
structural geometry (e.g., the span of the beam in the case of bending). 

Both explicit and implicit formulations of gradient plasticity have been investigated, with 
a second-order or fourth-order enhancement in the explicit case and with different types of 
boundary conditions in the implicit case. The gradient term incorporated into the softening 
law contains a new material parameter, which reflects the intrinsic material length scale. 

The results obtained with the second-order and fourth-order explicit model, as well as 
with the implicit model that applies homogeneous Neumann boundary conditions at the 
physical boundary of the body of interest, are qualitatively similar. Plastic yielding starts 
as soon as the yield stress is attained at the weakest section, the plastic zone continuously 
grows from that section and expands up to a maximum size that directly depends on the 
intrinsic material length and corresponds to the solution of the bifurcation problem for an 
idealized bar with perfectly uniform properties. At early stages of the inelastic process, the 
growth of the plastic zone occurs at increasing axial force (and thus at increasing stress). 
Even though the local material response is postulated as softening, the global response at the 
structural level is initially hardening, and only after a certain critical size of the plastic zone 
has developed, the global response turns into softening. The maximum level of the axial force 
thus depends not only on the yield stress and on the area of the weakest section but also on 
the distribution of the area in the vicinity of that section, which is in perfect agreement with 
the concept of nonlocal interactions in the material microstructure, described by the gradient 
terms. The increase of the ultimate load as compared to the local model (which would fail 
abruptly with plastic yielding fully localized into the weakest section) directly depends on 
the ratio between the length parameters characterizing the material and the geometry. If 
this ratio tends to zero, the ultimate load tends to the "locally" determined one. 

It is thus concluded that the behavior of the models mentioned in the previous paragraph 
is perfectly reasonable and all of them can describe a gradual development of the plastic zone, 
at the structural level accompanied by a transition from hardening to softening. In contrast 
to that, a modified version of the implicit gradient model with homogeneous Neumann 
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conditions prescribed at the boundary of the plastic zone (and not at the physical boundary 
of the body) exhibits a pathological behavior related to the unrealistic distribution of nonlocal 
plastic strain in the elastic zone, and it would lead to spurious yielding at points far from the 
"main" plastic zone, which are under low stress but their yield stress is artificially reduced. 
Therefore, this modified formulation of implicit gradient plasticity is not acceptable. 

In the future, the localization problem will be reformulated using a variational approach, 
which permits a more systematic treatment of discontinuities and reduces the regularity 
requirements. This is important e.g. for a rigorous justification of the admissibility conditions 
for the fourth-order gradient model, which are in this paper postulated based on intuition. 
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